This notebook is an exercise in the Geospatial Analysis course. You can reference the tutorial at this link.
Introduction
You are an urban safety planner in Japan, and you are analyzing which areas of Japan need extra earthquake reinforcement. Which areas are both high in population density and prone to earthquakes?

Before you get started, run the code cell below to set everything up.
import pandas as pd
import geopandas as gpd
import folium
from folium import Choropleth, Circle, Marker
from folium.plugins import HeatMapWe define a function embed_map() for displaying interactive maps. It accepts two arguments: the variable containing the map, and the name of the HTML file where the map will be saved.
This function ensures that the maps are visible in all web browsers.
Exercises
1) Do earthquakes coincide with plate boundaries?
Run the code cell below to create a DataFrame plate_boundaries that shows global plate boundaries. The “coordinates” column is a list of (latitude, longitude) locations along the boundaries.
plate_boundaries = gpd.read_file(r"C:\Users\LG PC\Desktop\data_mining\archive\Plate_Boundaries\Plate_Boundaries\Plate_Boundaries.shp")
plate_boundaries['coordinates'] = plate_boundaries.apply(lambda x: [(b,a) for (a,b) in list(x.geometry.coords)], axis='columns')
plate_boundaries.drop('geometry', axis=1, inplace=True)
plate_boundaries.head()| HAZ_PLATES | HAZ_PLAT_1 | HAZ_PLAT_2 | Shape_Leng | coordinates | |
|---|---|---|---|---|---|
| 0 | TRENCH | SERAM TROUGH (ACTIVE) | 6722 | 5.843467 | [(-5.444200361999947, 133.6808931800001), (-5…. |
| 1 | TRENCH | WETAR THRUST | 6722 | 1.829013 | [(-7.760600482999962, 125.47879802900002), (-7… |
| 2 | TRENCH | TRENCH WEST OF LUZON (MANILA TRENCH) NORTHERN … | 6621 | 6.743604 | [(19.817899819000047, 120.09999798800004), (19… |
| 3 | TRENCH | BONIN TRENCH | 9821 | 8.329381 | [(26.175899215000072, 143.20620700100005), (26… |
| 4 | TRENCH | NEW GUINEA TRENCH | 8001 | 11.998145 | [(0.41880004000006466, 132.8273013480001), (0…. |
Next, run the code cell below without changes to load the historical earthquake data into a DataFrame earthquakes.
# Load the data and print the first 5 rows
earthquakes = pd.read_csv(r"C:\Users\LG PC\Desktop\data_mining\archive\earthquakes1970-2014.csv", parse_dates=["DateTime"])
earthquakes.head()| DateTime | Latitude | Longitude | Depth | Magnitude | MagType | NbStations | Gap | Distance | RMS | Source | EventID | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 1970-01-04 17:00:40.200 | 24.139 | 102.503 | 31.0 | 7.5 | Ms | 90.0 | NaN | NaN | 0.0 | NEI | 1.970010e+09 |
| 1 | 1970-01-06 05:35:51.800 | -9.628 | 151.458 | 8.0 | 6.2 | Ms | 85.0 | NaN | NaN | 0.0 | NEI | 1.970011e+09 |
| 2 | 1970-01-08 17:12:39.100 | -34.741 | 178.568 | 179.0 | 6.1 | Mb | 59.0 | NaN | NaN | 0.0 | NEI | 1.970011e+09 |
| 3 | 1970-01-10 12:07:08.600 | 6.825 | 126.737 | 73.0 | 6.1 | Mb | 91.0 | NaN | NaN | 0.0 | NEI | 1.970011e+09 |
| 4 | 1970-01-16 08:05:39.000 | 60.280 | -152.660 | 85.0 | 6.0 | ML | 0.0 | NaN | NaN | NaN | AK | NaN |
The code cell below visualizes the plate boundaries on a map. Use all of the earthquake data to add a heatmap to the same map, to determine whether earthquakes coincide with plate boundaries.
# Create a base map with plate boundaries
m_1 = folium.Map(location=[35,136], tiles='cartodbpositron', zoom_start=5)
for i in range(len(plate_boundaries)):
folium.PolyLine(locations=plate_boundaries.coordinates.iloc[i], weight=2, color='black').add_to(m_1)
# Your code here: Add a heatmap to the map
HeatMap(data=earthquakes[['Latitude', 'Longitude']], radius=10).add_to(m_1)
m_1So, given the map above, do earthquakes coincide with plate boundaries?
2) Is there a relationship between earthquake depth and proximity to a plate boundary in Japan?
You recently read that the depth of earthquakes tells us important information about the structure of the earth. You’re interested to see if there are any intereresting global patterns, and you’d also like to understand how depth varies in Japan.
# Create a base map with plate boundaries
m_2 = folium.Map(location=[35,136], tiles='cartodbpositron', zoom_start=5)
for i in range(len(plate_boundaries)):
folium.PolyLine(locations=plate_boundaries.coordinates.iloc[i], weight=2, color='black').add_to(m_2)
# Your code here: Add a map to visualize earthquake depth
def color_producer(val):
if val <= 100:
return 'forestgreen'
else:
return 'darkred'
# Add a bubble map to the base map
for i in range(0,len(earthquakes)):
Circle(
location=[earthquakes.iloc[i]['Latitude'], earthquakes.iloc[i]['Longitude']],
radius=20,
color=color_producer(earthquakes.iloc[i]['Depth'])).add_to(m_2)
# Display the map
m_2Can you detect a relationship between proximity to a plate boundary and earthquake depth? Does this pattern hold globally? In Japan?
3) Which prefectures have high population density?
Run the next code cell (without changes) to create a GeoDataFrame prefectures that contains the geographical boundaries of Japanese prefectures.
# GeoDataFrame with prefecture boundaries
prefectures = gpd.read_file(r"C:\Users\LG PC\Desktop\data_mining\archive\japan-prefecture-boundaries\japan-prefecture-boundaries\japan-prefecture-boundaries.shp")
prefectures.set_index('prefecture', inplace=True)
prefectures.head()| geometry | |
|---|---|
| prefecture | |
| Aichi | MULTIPOLYGON (((137.09523 34.65330, 137.09546 … |
| Akita | MULTIPOLYGON (((139.55725 39.20330, 139.55765 … |
| Aomori | MULTIPOLYGON (((141.39860 40.92472, 141.39806 … |
| Chiba | MULTIPOLYGON (((139.82488 34.98967, 139.82434 … |
| Ehime | MULTIPOLYGON (((132.55859 32.91224, 132.55904 … |
The next code cell creates a DataFrame stats containing the population, area (in square kilometers), and population density (per square kilometer) for each Japanese prefecture. Run the code cell without changes.
# DataFrame containing population of each prefecture
population = pd.read_csv(r"C:\Users\LG PC\Desktop\data_mining\archive\japan-prefecture-population.csv")
population.set_index('prefecture', inplace=True)
# Calculate area (in square kilometers) of each prefecture
area_sqkm = pd.Series(prefectures.geometry.to_crs(epsg=32654).area / 10**6, name='area_sqkm')
stats = population.join(area_sqkm)
# Add density (per square kilometer) of each prefecture
stats['density'] = stats["population"] / stats["area_sqkm"]
stats.head()| population | area_sqkm | density | |
|---|---|---|---|
| prefecture | |||
| Tokyo | 12868000 | 1800.614782 | 7146.448049 |
| Kanagawa | 8943000 | 2383.038975 | 3752.771186 |
| Osaka | 8801000 | 1923.151529 | 4576.342460 |
| Aichi | 7418000 | 5164.400005 | 1436.372085 |
| Saitama | 7130000 | 3794.036890 | 1879.264806 |
Use the next code cell to create a choropleth map to visualize population density.
# Create a base map
m_3 = folium.Map(location=[35,136], tiles='cartodbpositron', zoom_start=5)
# Your code here: create a choropleth map to visualize population density
Choropleth(geo_data=prefectures['geometry'],
data=stats['density'],
key_on="feature.id",
fill_color='YlGnBu'
).add_to(m_3)
# Display the map
m_3